% funcion [If]=reinicializar(I, max_iter , delta_t, tol, [verbose], [tipo]),
% realiza iteraciones con paso delta_t de la discretización de la
% ecuacion de reinicializacion: u1=u0-delta_t*(sign(I).*(sqrt(u0_x.^2+u0_y.^2)-1))
% Input:
% I: superficie original, (con la linea de nivel 0 asociada a la curva)
% max_iter: número de iteraciones maximo
% delta_t: paso
% tol: toleracion
% verbose:[default 0]
% tipo: tipo de discretizacion
%       - 
%       - 
% Output:
% If: funcion distancia que mantiene la linea de nivel 0 que la original

function [If]=reinicializar(I, max_iter, delta_t, tol, varargin)

%% inicializacion de las variables y valores por defecto. --------
if nargin<5,
    verbose=0;
else
    verbose=varargin{1};
end

if nargin<6,
    tipo='def';
else
    tipo=varargin{2};
end
% ---------------------------------------------------------------

%% inicializacion
[m,n]=size(I);
u0=I;
dif=2*tol;
mostrar_cada=30;
num_iter=0;
tipo_gradiente='centrada'; con_bordes='nulas';

if verbose, %abro las ventanas que necesito
    h1=figure;
end

% ---------------------------------------------------------------

%% loop,
while num_iter<max_iter && dif>tol,
    num_iter=num_iter+1;
    
    switch tipo,
        case 'def',
            aux=grad(u0,con_bordes,tipo_gradiente);
            u0_x=aux(:,:,1); u0_y=aux(:,:,2);
            u1 = u0 - delta_t * ( sign(I) .* (sqrt(u0_x.^2+u0_y.^2)-1) );
        case 'paper'
            ap=zeros(m,n); bp=ap; cp=ap; dp=ap;
            an=zeros(m,n); bn=an; cn=an; dn=an;
            
            % variables auxiliares
            ap(2:end-1,:)=max(u0(2:end-1,:)-u0(1:end-2,:),0);
            bp(2:end-1,:)=max(u0(3:end,:)-u0(2:end-1,:),0);
            cp(:,2:end-1)=max(u0(:,2:end-1)-u0(:,1:end-2),0);
            dp(:,2:end-1)=max(u0(:,3:end)-u0(:,2:end-1),0);
            
            an(2:end-1,:)=-min(u0(2:end-1,:)-u0(1:end-2,:),0);
            bn(2:end-1,:)=-min(u0(3:end,:)-u0(2:end-1,:),0);
            cn(:,2:end-1)=-min(u0(:,2:end-1)-u0(:,1:end-2),0);
            dn(:,2:end-1)=-min(u0(:,3:end)-u0(:,2:end-1),0);
            
            epsilon=1;
            sgnI=I./sqrt(I.^2+epsilon^2);
            
            G_u0=zeros(m,n);
            
            G_u0=(I>0) .* ( sqrt( max(ap.^2 , bn.^2) + max(cp.^2 , dn.^2) ) - 1 ) +...
                (I<0) .* ( sqrt( max(an.^2 , bp.^2) + max(cn.^2 , dp.^2) ) - 1 );
            
            u1 = u0 - delta_t * (sgnI .* G_u0);
                        
        otherwise
            disp('ERROR: opcion en "reinicializar" invalida')            
    end
    
    dif = mean( abs( u1(:)-u0(:) ) );
    
    u0 = u1;
    
    if verbose,
        disp(['[Ec Reinicializacion] Iteracion número ' num2str(num_iter) ' diferencia = ' num2str(dif) ])
        if mod(num_iter,mostrar_cada)==0,
            figure(h1),
            set(h1,'name',['u0 en el tiempo: ' num2str(num_iter*delta_t)],...
                'NumberTitle','off');
            surfc(u0),hold on, surf(0*u0), hold off,...
            shading interp, view(120,2),
 
        end
    end
    
end


if num_iter==max_iter,
    %disp('OJO: se alcanzo el maximo de iteraciones en la pde auxiliar.')
end

If=u0;

%% Cierro las que abri,
if verbose, 
    close(h1),
end
        
    